Collective Oscillations in Classical Nonlinear Response of a Chaotic System 
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We consider classical response in a strongly chaotic (mixing) system. As opposed to the case of 
stable dynamics, the nonlinear classical response in a chaotic system vanishes at large times. The 
physical behavior of the nonlinear response is attributed to the exponential time dependence of the 
stability matrix. The response also reveals certain features of collective resonances which do not 
correspond to any periodic classical trajectories. We calculate analytically linear and second-order 
response in a simple chaotic system and argue on the relevance of the model for interpretation of 
spectroscopic data. 
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Time-domain femtosecond spectroscopy constitutes a 
powerful tool that probes electronic and vibrational co- 
herent dynamics of complex molecular systems in con- 
densed phase ^.i^:^iii^. Spectroscopic signals are directly 
related to optical response functions that carry detailed 
information on the underlying dynamical phenomena. At 
room temperatures the complexity often originates from 
slow strongly anharmonic vibrational modes that can be 
treated within the framework of classical mechanics. 

A number of studies have been devoted to the classi- 
cal response in stable (integrable) dynamical systems^i^ii^. 
The nonlinear response functions have been shown di- 
verge linearly with time, while the divergence can be 
eliminated by invoking a fully quantum description^'S'. 
Although integrable dynamics represents only an approx- 
imation for realistic physical situations, weak deviations 
from the integrability do not eliminate the unphysical be- 
havior of the response functions. This follows from the 
fact that the divergence is related to quasiperiodic mo- 
tions on invariant tori§. According to the Kolniogorov- 
Arnold-Moser (KAM) theory, most the invariant tori are 
not destroyed by small perturbations that break down 
integrability^. Yet, stable dynamics would be typical 
for a close to equilibrium situation. At larger energies a 
generic situation would correspond to dynamical behav- 
ior with chaotic features^. Moreover, unstable (hyper- 
bolic) dynamics may be more common due to the stabil- 
ity of chaos with respect to perturbations. It has been ar- 
gued based on results of numerical analysis-!^ that chaotic 
dynamics appears to observe the convergence of the clas- 
sical response functions. In spite of apparent importance, 
and to the best of our knowledge the problem of the non- 
linear response in strongly chaotic systems has never been 
addressed using analytical methods. We note that ana- 
lytical calculations of the response are rarely feasible in 
nonintegrable systems, whereas numerical simulations of 
chaotic dynamics are complicated by the exponential di- 
vergence of stability matricesiS. 

In this Letter we show that (i) the classical response 
of the chaotic system exhibits decay and oscillations as a 
function of times between probing pulses, and (ii) the 



Fourier transform of 2D second order response func- 
tion reveals broad and asymmetric peaks as signatures 
of chaos. This is the main result or our study. 

A strongly chaotic (mixing) system is characterized by 
a special spectrum of the Liouville operator^''' The 
spectrum consists of complex Ruelle-Pollicott (RP) reso- 
nances that determine the asymptotic oscillations and 
decay of the correlations. This yields the linear re- 
sponse function directly related to the two-point correla- 
tion functions, in agreement with the fluctuation dissipa- 
tion theorem (FDT). Nonlinear response, however, turns 
out to be more involved still with the noticeable effect of 
the resonances. 

A system driven by a time-dependent external field 
£{t) is described by the Hamiltonian Ht = H — f£{t) 
with the the function /(rj) in phase space representing 
the dipole (see, e.g. Ref. 111). The response functions 
S'(") that depend on n time intervals describe the ex- 
pansion of the measured signal (/(i7(i))) — J dt] pf in a 
functional series in £. The second-order response func- 
tion reads 

S^^\hM)^dt,j dr,fe~^'-[f,e-^'-fdEPo]. (1) 

where po is the equilibrium distribution and L = {H, •} 
is the Liouville operator of the unperturbed system. 

A schematic picture of second-order response forma- 
tion in a chaotic system that employs the Liouville phase 
space evolution is shown in Fig. [1] Propagating the 
observable / in Eq. ([1]) backward in time we interpret 
the second-order response as an overlap of the distri- 
butions f-{ri) — eyip{Lt2)f and £^^djf+{r]) with /+ — 
exp{—Lti)fdEPo and the vector field S^^dj = {/, •}. Since 
/ dxf = 0, we can for simplicity represent the function 
/ by two separate regions in the phase space where it 
adopts positive and negative values. In the case of hy- 
perbolic dynamics evolution during time t changes the 
regions' shapes. For \t\ ^ 1, the shape becomes sim- 
ilar to ribbon-like fettuccine: elongated along the un- 
stable direction by a factor ^ e'^*, narrowed along the 
stable one ~ e~^* (with A being the Lyapunov exponent) 
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FIG. 1: Schematic cross-section of the phase space along the 
surface given by the stable and unstable directions. Initial 
distribution of / is presented by two regions + and — (dark 
red and light green). As time elapses the distribution elon- 
gates along unstable direction and contracts along stable one. 



and unchanged along the flow. Since /_ results from 
reverse dynamics, the distributions /+ and /_ are elon- 
gated along the unstable and stable directions, respec- 
tively. Therefore, the overlap is represented by a large 
set of iV ~ gHti+t2) small disconnected regions with the 
volume V ~ e"''''^*^^*^''. Since / is represented by a posi- 
tive and negative regions, the distribution /_ consists of 
two positively and negatively "charged" fettuccine and 
the cancellations result in a signal determined by a typi- 
cal fluctuation proportional to \/]V, and the overlap inte- 
gral attains a factor y/Nv ^ g-Hti+t2)/2 -j-j^g^-j. ^^^ns out 
to be exponentially small. Yet, this is not the end of the 
story since the derivative S,^dj of a sharp feature along 
the stable direction can create exponentially large ~ C^*^ 
factors. This is the Liouville space signature of the ex- 
ponentially growing components of the stability matrix, 
which affects the response starting with second order due 
to FDT^'^'^^. However, the divergent terms cancel out: 
decomposing ^ = + ^-i- + ^- into the direction along 
the flow, and unstable and stable components, we note 
that only the last term is potentially dangerous. Cal- 
culating the dangerous component of the overlap inte- 
gral by parts we arrive at two contributions: the over- 
lap of ^■'_djf- with /+, and the overlap of /_ with /+ 
weighted with div^_ . The first contribution contains an 
additional e~^*'^ factor, the second one provides an ad- 
ditional factor independent of time. This results in a 
physical ~ exp~''^(*i+*^^/^ large time asymptotics of the 
nonlinear response function. 

Trajectories of a particle in arbitrary potential U{r) 
are known to be the same as for a free motion in a curved 
space with the metric gik = (1 — U{r) / E)5ii~^. Al- 
though the potential generates nonuniform motion along 
the trajectories, one can expect chaotic behavior due to 
the exponentially growing separation between the trajec- 
tories. In addition, when the motion is finite, the acces- 
sible part of the configurational space at a given energy 



can be multiply connected. In the simplest case of two 
coordinates the motion occurs inside a disk-like region 
punctured by g forbidden islands. Some fraction of tra- 
jectories approaches the boundaries so close that this can 
be qualified as reflection. Utilizing the original argument 
of Sinaii^ reflection can be interpreted as continuing mo- 
tion on the antipode replica of the accessible region glued 
to its original coun terp art via the boundaries (see Ap- 
pendix [X] and Ref. [iO). The resulting compact surface 
has topology of a sphere with g handles (Riemann surface 
of genus g). In the 5 > 1 case the average curvature is 
negative, which results in unstable (hyperbolic) dynam- 
ics. 

To rationalize the described above qualitative picture 
of response in a chaotic system we calculate the linear and 
second-order classical response functions for free motion 
in a Riemann surface of constant negative (Gaussian) 
curvature. This model that allows an exact solution has 
been serving as a prototype of classical chaos, as well as 
an example of semiclassical quantization^ i ^° ' ^ ^ . 

The classical free-particle Hamiltonian 

H{x, C) = (l/2m)g'V<Pfe = CV2m , (2) 

depends on the absolute momentum value C, only, and 
X G includes the two coordinates and momentum 
direction angle 9. Therefore, the smooth compact 3-D 
manifold , which represents the subspace of phase 
space with fixed energy, is preserved in classical dynam- 
ics. Hereafter, we will use dimensionless units so that 
m = 1 the curvature K = —1, and C = 1 (when energy 
is fixed). Our model allows for an exact solution due 
to strong Dynamical Symmetry (DS). To describe DS we 
adopt an agreement used in differential geometry by iden- 
tifying a first-order differential operator of differentiating 
along the vector field with the vector field itself. We de- 
note by (Ti the vector field that describes the phase space 
dynamics, and set Uz — d/d9, ai = [cri, Uz]- A simple lo- 
cal calculation yields [1J2, CTz] = —ai and [cti, 02] = —Ka^ 
which implies that in the constant negative curvature 
case the vector fields cr^, ui, and = [cri,az] form the 
Lie algebra so(2, 1). The group 50(2, 1) action in the re- 
duced phase space is obtained by integrating of the 
so(2, 1) algebra action. DS with respect to the action 
of the group G = SO (2, 1) does not mean symmetry in 
a usual sense, i.e. that the system dynamics commutes 
with the group action, but rather reflects the fact that 
the Poisson bracket is given by a bi-vector field 

Lu = d( iSi ai - <Ti iSi d( + ((T2 <Siaz ~ (Tz® 0-2) . (3) 

In particular, the vector field that determines classical 
dynamics is given by an element cti of the correspond- 
ing Lie algebra so(2, 1), whereas the stable and unstable 
directions are determined by T "'2 • 

DS implies that the space of phase-space distributions 
constitutes a representation of SO{2, 1) and, being de- 
composed into a sum of irreducible representations, pro- 
vides a set of uncoupled evolutions. Unitary irreducible 
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representation of 5*0(2, 1) are well-known and can be 
conveniently implemented in terms of functions ^(u) in 
a circle^. Principal series representations relevant for 
our calculations are labeled by imaginary parameter s 
with the Liouville operator cti = sinu9„-|-(l/2 — s) cosu. 
The aforementioned decomposition identifies the angu- 
lar harmonic ^'fe(u) = e*'^" in a circle with a phase- 
space distribution ijjkix^s) with given angular momen- 
tum azipk{x;s) = iki]jk{x]s). Commutation law in 
so(2, 1) entails <j±'il:k{x\ s) = (±fc -I- 1/2 — s)'il}k±i{x] s) 
with two anti-Hermitian conjugated ladder operators 
a± — ai zL ia2, ct^± — Tip- Any relevant distribution can 
be decomposed in the modes ipk{x;s), where s adopts 
discrete values (the spectrum). The spectrum {sq} can 
be related to the eigenvalues of the Laplacian opera- 
tor in by identifying = — ^((t+(t_ -I- cr_cr+) and 

-V2^o(a;;s) = (l/4-s2)^o(a;;s). 

The dipole / is a function in Af^, i.e. a function in 
phase space that docs not depend on both and mo- 
mentum direction 9, hence a^f — 0, and / can be ex- 
panded as a sum over the principal series representations 
/ = J2sBs^oix;s). 

The integrands in the response function ([1]) include 

the result of action of the evolution operator e^^* on the 
components ipo{x;s) of the dipole momentum and their 
derivatives resulting from the Poisson brackets. To sim- 
plify the calculations, we transform the expressions of 
response functions in a way that the evolution of deriva- 
tives of ipo{x', s) does not appear. 

We use the representation in the circle to find the ex- 
pansion of e~^*ij}Q{x; s) over basis vectors ipk{x] s), 
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Akit] s)ipk{x; s) 



(4) 



Since e~^*'I'o(M) represents a function obtained as 
the result of the action of the evolution operator on 
^o{u) = 1, it can be found by solving the equation 
dtg{t,u) + Lg{t,u) = supplemented with initial con- 
dition g(0,u) = 1. In the principal series representation 
with Res = 0, Ims > 0, and the equation takes the 
following form: 

dtg{t, u) + {sinudu + (1/2 — s) cosm) g{t, u) = . (5) 

The solution can be easily found, and thus after the in- 
tegration J due~^^^g{t,u) we obtain Ak{t]s): 



Ak[t\s) 



Re 



r(s)e^ 



^/^^(l/2 - s) 
1 



T{k+\ + s) 



iFi k 



2 -s,A;+-,l-s,e" 



where 2F1 is the Gauss hypergeometric function. For 
details see Appendix IB] and Ref. 21. 

The linear response function can be obviously ex- 
pressed via A[){t]s). The calculation involves only the 



first part of the Poisson bracket, and the result is 
conveniently represented similarly to the fluctuation- 
dissipation relation S'^^\t) cx dtA^lt] s)^^. 

For large t the linear response function shows damped 
oscillations e^^'*"^/^-'*. The expansion in powers of e~^* 
corresponds to RP resonances^'^. Only even RP reso- 
nances contribute to the response function. The expan- 
sion is a converging series in e~^* if e~^* < 1, i.e. i > 0. 
Since the dipole f{x) can decomposed in irreducible rep- 
resentations characterized by different values of s, the 
linear response function constitutes a linear combination 
of contributions with coefficients B^. 

The second-order response function can be calculated 
using Eqs. ^ and ([3]^. The calculation can be 
substantially simplified by propagating the observable 
/ in Eq. ([T]) backwards in time and making use of 

^~Lt2^ _ e^*^. Then we apply Eq. ^ to decompose 

e-^*^ Vo(3;; s) in %ljk{x\s) and e^-^*^ Vo(a^; in ^]i{x;s). 
The integration over the reduced phase space includes 
an integral over the momentum direction 9 that results 
in vanishing of all terms with k ^ / ± 1. The second 
order response function consists of several contributions, 
5^'^ = E BpBqBrY.V=Q{~'^fsfl of similar form, 
e.g.: 



with the matrix elements 



dx ipl{x; si)V'o(a;; S2)^pk{x\ S3) 



(6) 



We further establish recurrent relations for coefficients 
([S]) that allow expressing all of them via few "initial 

For instance, in the sim- 
s, the recurrent relations 



e.g. 



conditions" 
plest case si = S2 = S3 
read Ofc+i = TUxfef^afc 



(2fc-|-l)^-s^"fc (2fe+l)^ „ 

efficients are expressed via oq due to the symmetry 
ak — a-k^^- In this case we find the asymptotic behav- 
ior ak oc k^^/^^^ for k — s- 00. The summation over k 
is performed numerically. The set of coefficients ao, as 
well as the spectrum {Sq,} are attributes of a particular 
Riemann surface. Constant curvature Riemann surfaces 
of genus g are classified by the so-called moduli spaces 
whose dimensions grow linearly in g. Therefore, any fi- 
nite set Sa of spectral elements, as well as the relevant 
set of coefficients ag can be implemented for some partic- 
ular Riemann surface, and hereafter we will treat them 
as independent parameters. 

The convergence of the series over k is ensured by the 
dependence Ak{t; s) oc exp(— 2A:e~*) for large k. The se- 
ries is almost sign-alternating, and the dependence of 
the summand magnitude on k becomes smoother with 
increasing k, independently of t. This allows for an ef- 
fective procedure that evaluates the series, since it does 
not involve numerical summation of the exponentially in- 
creasing with t number of terms. The procedure also 
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works in any order of e~*^ and e~*^ obtained from the 
hypergeometric expansions of Af^{t; s)^. In Ref. ^ we 
extend our study of optical response by adding Langevin 
noise to classical deterministic dynamics. We identify the 
RP resonances as the eigenvalues of the related Fokker- 
Planck operator C ~ — kV^ + L in the weak-noise limit 
K ^ 0. We demonstrate that a standard spectral decom- 
position of the response functions in the eigenmodes of 
the Fokker-Planck operator C has a well-defined limit at 
K —^ that reproduces our expression for the classical 
response functions. 
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FIG. 2: Absolute value of 2D Fourier transform of the second 
order response function: (a) single resonance s — 5i, (b) linear 
combination of terms with two resonances si = 3* and S2 — 
5i. Linear plots show cross-sections of the spectra at uji — 
UJ2 = 5. 

The absolute value of the 2D Fourier transform of the 



second-order response function presented in Fig. [2] shows 
diagonal and cross peaks, as well as a stretched along 
the LUi direction feature (real and imaginary parts are 
presented in Appendix [Cjl . The latter originates from 
time-domain damped oscillations with variable period 
and can be interpreted as a signature of chaos (insta- 
bility) in the underlying dynamics. Also note that in our 
chaotic case the peak frequencies may not be attributed 
to any particular periodic motions, although they can be 
expressed in terms of all periodic orbits via the dynam- 
ical C- function^"* . Therefore, they can be referred to as 
collective chaotic resonances. 

Our results can be applied to the interpretation of spec- 
troscopic data. Classical chaos is quite generic for dy- 
namics in large molecules and molecular networks cou- 
pled via hydrogen bonds. Even for a small number of 
vibrational modes, e.g. in small systems of hydrogen 
bonds^"^, the shape of the effective potential energy can 
lead to chaotic dynamics. Moreover, the chaos may be 
caused by just some regions with nonuniform curvature 
of the potential surface^i. 

In the present Letter we studied classical nonlinear re- 
sponse in a strongly chaotic (mixing) system. Detailed 
analysis performed for free motion in a Riemann com- 
pact surface of constant negative curvature demonstrated 
that calculation of the response based on Liouville space 
dynamics does not have fundamental difficulties. We 
demonstrated that peaks in 2D spectra are not necessar- 
ily attributed to periodic or quasiperiodic motions, but 
rather can have collective nature. We argue that col- 
lective resonances should be accompanied by stretched 
features in 2D spectra that constitute a spectroscopic 
signature of underlying dynamical chaos. 
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APPENDIX A: EFFECTIVE NEGATIVE 
CURVATURE OF CONFIGURATION SPACE 

In this appendix we rationalize a close qualitative 
connection between classical dynamics of multidimen- 
sional motion in a potential with forbidden islands and 
a geodesic flow in a compact manifold with non-trivial 
topology. For the sake of simplicity we consider the case 
of two coordinates. The motion is restricted to a disk-like 
region punctured by forbidden islands. Boundaries of the 
accessible region correspond to the lines where the total 
energy coincides with the potential energy. To employ 
the original argument of Sinaii^ we glue the accessible 
region of the configuration space to its antipode replica 
along the boundary. A reflection from the boundary can 
be thought of as continuing motion in the other compo- 
nent. The resulting surface is compact and its topology is 
characterized by genus g (the number of handles attached 
to a sphere) corresponding to the number of original is- 
lands. In Fig. [3] we show both the original billiard-like 




FIG. 3: Motion in the multiply connected flat configura- 
tion space is equivalent to motion on the curved surface. 
Idealization of the hard-wall potential: Trajectories are con- 
fined within the original classically accessible 27? configura- 
tion space reflecting on its boundaries (left). Reflections can 
be viewed as transitions to the antipode surface component 
glued to the original one along the boundaries (right). Trajec- 
tories on the deformed smooth version of the resulting com- 
pact surface are shown as solid and dashed lines, which corre- 
spond to the original and antipode components, respectively. 



motion in the configuration space and the equivalent mo- 
tion in two components assembled together that form a 
compact surface of genus g = 2. As mentioned in the 
main text trajectories of a particle moving in a poten- 
tial resemble the geodesic lines in curved space with the 
metric gik = (1 — U[r)/ E)5ik- The regions of the con- 
figuration space that correspond to negative curvature 
work as a defocusing lenses, causing instability (diver- 
gence of close trajectories). For example, close trajecto- 
ries diverge every time they approach the boundary of the 
classically inaccessible island or pass a region of a poten- 
tial local maximum that belongs to the accessible region. 
Passing through the stable regions with positive curva- 
ture cannot compensate for the instability, since stabil- 
ity reflects oscillatory, rather than converging features in 
the dynamics of the trajectory deviation. In particular, 
existence of unstable regions combined with ergodicity 
ensures exponential divergence of trajectories over long 
enough times. According to the Gauss-Bonnet theorem, 
in the g > 1 case the average curvature is negative, which 



implies regions of instability. 



APPENDIX B: LIOUVILLE-SPACE 
CALCULATION OF RESPONSE FUNCTIONS 

Response functions can be obtained starting with the 
driven Liouville equation dtp+C = £{t){f, p}. The phase 
space distribution p(r], t) is represented as a functional se- 
ries in £{t) with the initial distribution po being the zero- 
order contribution. In our case of interest the observable 
coincides with the dipole moment / that represents cou- 
pling to the probe field. The linear and second-order 
classical response functions are expressed via the Poisson 
brackets as 

^W(i)- J dvf{v)e-^'{fiv),Po}, (Bl) 

S^'HhM) = I dvfiri)e-^'Hf{v),e~^'Hfiv),Po}}- 

We start with utilizing the fact that the equilibrium dis- 
tribution Po depends on energy E = (^'^ /2 only, and make 
a transformation based on the FDT. This brings us to 
Eq. (P). 

The problem of free motion in a compact sur- 
face of constant negative curvature possesses dynamical 
SO {2, 1) symmetry. DS means that the flow is deter- 
mined by the so(2, 1) generator cti, it also conserves the 
stable and unstable directions; this is expressed by the 
relation e'^^\a2 ± cr^) = e^*(cr2 ± (j^)e'^i*. 

Any phase-space distribution can be decomposed into 
a direct sum of irreducible representations; only principal 
series representations of 50(2, 1) contribute to the linear 
and second-order response^^-. DS implies that the distri- 
butions in different representations evolve independently. 
Distributions in irreducible representation can be imple- 
mented as functions in a circle with the basis of angular 
harmonics '^k{u) = e'*^". The evolution in the circle is 
determined by Eq. ^ that, given the initial condition 
g{0,u) = \E'(u) = 1, can be solved using the method of 
characteristics: 

g(t, u) = (cosh i -f sinh< cos u)" . 

Implementing the correspondence of natural scalar prod- 
ucts in and in the circle, we find the expansion coef- 
ficients in Eq. (|4]) 

Akit;s)^ J ^e-^'"git,u). 



These coefficients represent the main dynamical ingredi- 
ent in the response function calculation. We represent 
them in terms of the hypergeometric functions. Alterna- 
tive representations have been derived in the context of 
two-point correlationsiiii^. The linear response function 
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is determined by the dynamical part alone: 



oo 



Calculation of the second-order response function is 
more involved. Substituting the Poisson brackets ^ into 
Eq. (|B1|) the 2D signal adopts a form 

+ (e-^*V)*(^2/)(a.e--^«*^/)}^. 

All terms in parentheses are further expanded in irre- 
ducible representations and angular harmonics. The re- 
sulting contributions include two coefficients A^. and a 
triple product given by Eq. ©. The latter is the ge- 
ometrical ingredient that appears in nonlinear response 
functions. The simplest expression for the second-order 
response function that incorporates different collective 
resonances, originates from two terms in the decomposi- 
tion / = Bs-^ipo{x; si) + Bg^ipoi^] 82)- Then the response 
function reads 



p,q,r—si ,S2 



ri=0 



n + - + r \ a- 



■ d 

- n 



{Al{t2;p)An+i{ti;r)) 



''n+l 



Ot2 



dynamical system considered in this Letter. Experimen- 
tal data on 2D time-domain spectroscopy that probes 
the response function S^'^\ti,t2) is usually presented us- 
ing the so-called 2D spectra that constitute a numerical 
2D Fourier transform of the response function with re- 
spect to ti and ^2- The absolute value of the 2D spectra 
for our dynamical system are given in Fig. [21 whereas in 
Fig. 2] we present the real and imaginary parts. 
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FIG. 4: Real and imaginary parts of 2D spectra: (a) real 
part, (b) imaginary part with the single resonance s — 5i; (c) 
real part, (d) imaginary part of linear combination of terms 
with two different resonances si = 3i and S2 = 5i. 



It turns out that the resulting second-order response 
function may be represented as a double expansion over 
RP resonances. The expansion is nontrivial: one encoun- 
ters diverging series if the expansion of quantities An(t; s) 
is performed first. Fortunately, we have found an efficient 
scheme for numerical evaluation of the series^. The pro- 
cedure also demonstrates that the second-order response 
is indeed expanded in the RP resonances. 



APPENDIX C: 2D SPECTRA FOR 
SPECTROSCOPY OF A CHAOTIC SYSTEM 

In this appendix we present some details on 2D spectra 
that characterize the second-order response of a chaotic 



As follows from the general picture presented in the 
main text, for a given unperturbed dynamics the re- 
sponse function depends on the expansion of the dipole 
function / in the eigcnmodes of the Laplace operator: 
each eigenmode provides with an oscillation in the sig- 
nal. Since the dipole is generally a smooth function of 
the system coordinates, the number of eigenmodes that 
participate in the expansion is typically small. To study 
the important features of the signals we consider the cases 
of one and two oscillation in the signal. In the first case 
we see a diagonal peak with a stretched feature along wi 
direction; such stretching never occurs due to a periodic 
orbit in an integrable system even if it is coupled to a 
harmonic bath. In the second case we see both diagonal 
and cross peaks accompanied by the stretched features. 



